Segregation induced by inelasticity in a vibrofiuidized granular mixture 
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We investigate the segregation of a dense binary mixture of granular particles that only differ in 
their restitution coefficient. The mixture is vertically vibrated in the presence of gravity. We find a 
partial segregation of the species, where most dissipative particles submerge in the less dissipative 
ones. The segregation occurs even if one type of the particles is elastic. In order to have a complete 
description of the system, we study the structure of the fiuid at microscopic scale (few particle 
diameters). The density and temperature pair distribution functions show strong enhancements 
respect the equilibrium ones at the same density. In particular, there is an increase in the probability 
that the more inelastic particles group together in pairs (microsegregation). Microscopically the 
segregation is buoyancy driven, by the appearance of a dense and cold region around the more 
inelastic particles. 

PACS numbers; 45.70.-n, 45.70.Mg 



I. INTRODUCTION 

Granular matter, when externally excited, show a se- 
ries of peculiar phenomena. One of them is the mixing 
or segregation that takes place when two or more species 
of different grains are put together. Depending on the 
control parameters and the energy injection mechanisms, 
grains of different size, shape, mass or mechanical prop- 
erties can mix or segregate. 

Consider the particular case of a mixture of small 
grains and one large intruder. Here the intruder can go 
up [l| or down 3 -the direct and reverse Brazil nut ef- 
fects, respectively-, phenomenon that has been studied 
in many papers (see, e.g. Ref. Q and references therein). 
When both species have similar sizes (but possibly differ- 
ent) we can select few cases where a variety of segregation 
mechanisms and scenario appear For instance, in 

Ref. 0] particles of different masses, radii and restitution 
coefficients are placed in a dish which is horizontally vi- 
brated, finding complete segregation. Segregation is also 
found in the same geometry when the grains have differ- 
ent friction coefficient with the base fsl . Under horizontal 
swirling, radial segregation of particles of different sizes 
has been observed In avalanches, grains of different 
shape segregate in stripes jioj ; in partially filled rotating 
drums, axial size segregation develops In two di- 

mensional systems under gravity, sinusoidally vibrated, 
clustering has been observed ■ This segregation effect 
can be modulated by using non-sinusoidal vibration [l^ . 

In some of the cases mentioned above the grain species 
differ on the friction or restitution coefficient. However 
few papers have studied segregation when this is the only 
difference between grains. One of these cases is Ref. 
[T^ , where a mixture of spheres that only differ in fric- 
tion coefficients (static, dynamic and rolling) is horizon- 
tally vibrated. They find complete mixing -that is, no 
segregation- for a flat plate while segregation is only ob- 
served when the plate was slightly inclined. Therefore, 
these results contradict the previously mentioned ones. 



In a theoretical approach, Ref. constructs the 

hydrodynamic equations from the Boltzmann equation, 
finding segregation induced by inelasticity. The authors 
explain the phenomenon as a consequence of the tem- 
perature gradient in the system induced by inelastic col- 
lisions, and relate the concentration gradient with the 
temperature gradient. In the same spirit, Ref. [lB| stud- 
ies the low density hydrodynamics of a mixture in the so 
called tracer limit, i.e. where the concentration of one of 
the components tends to zero. Among other results, they 
find that the temperature ratio of both species must be a 
constant. This constant value was already measured by 
two experimental groups [3 in 2D and 3D respec- 
tively and by computer simulations (Toj . Generalization 
to high density has been done by Garzo f20j . 

The goal of this paper is to confirm or deny the exis- 
tence of segregation induced by a inelasticity difference 
and characterize this phenomenon. The main tool will 
be Molecular Dynamics computer simulations of two di- 
mensional systems of a binary mixture kept fluidized by 
a vibrating base. 

The structure of this paper is as follows. In Sect. II we 
describe the system under consideration. Section III is 
devoted to the macroscopic study of the system, in partic- 
ular density and temperature profiles. Section IV present 
a microscopic study via the pair distribution functions. 
Section V proposes a model that possibly explains the 
segregation. We conclude with Section VI summarizing 
the results of the paper. 



II. DESCRIPTION OF THE SYSTEM 

We study the effect of the difference on restitution co- 
efficients in the segregation phenomenon, by means of 
Molecular Dynamics simulations of a bidimensional gran- 
ular mixture of two types of particles, named A and B. 
Grains are modeled as smooth Inelastic Hard Disks both 
having the same diameter a and mass m, but differing on 
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the normal restitution coefficient that characterizes their 
inelastic collisions. The restitution coefficient for A-A 
collisions is a a, for B-B collisions is a^. For the inter- 
particle collisions A-B we have taken oab — (o^A+ctB)/'^- 
In what follows we will consider that B are the most in- 
elastic particles (as < ua)- 

We have taken a fixed total number of particles N-r = 
Na + Nb, changing the concentration of the B particles. 
For the simulations reported in this paper, we have fixed 
Nt = 680 disks and varied Nb from 10 (that can be 
considered as a tracer limit) until 160. The disks are 
placed under the action of a gravitational acceleration g 
pointing downward in a rectangular box of width = 
50(7, infinite height, and with the bottom wall oscillating 
periodically at high frequency to and small amplitude A, 
with a bi-parabolic waveform [2]| . Periodic boundary 
conditions are used in the horizontal direction, trying to 
avoid the appearance of convective rolls by the influence 
of the walls. Under these conditions, the system reaches 
a stationary state with gradients in the vertical direction 

Units are chosen such that cr = 1, m = 1, and we 
fix the energy scale by the wall oscillation, m{A(jj)'^ = 1. 
Simulations are performed with g = 0.15, and A = O.Olcr. 



III. MACROSCOPIC SEGREGATION 

In order to illustrate the main observed features, we 
report results of a simulation having a small fraction of 
inelastic particles Nb = 10, and = 0.7 and the rest 
nearly elastic: Na — 670, aA — 0.98. 

The density profiles of the two species are shown in Fig. 
[T^, where we plot the number density of particles of type 
A and B: r7,^(z) and nB{z). The normalization of these 
quantities is such that dz ^ dxnA{z) = Na (resp. 
B). For plotting purposes only, ub is rescaled by a fac- 
tor Na/Nb, so in the case of no segregation both proffles 
would be identical. Both densities have the characteristic 
shape of vibrofluidized systems subject to gravity: there 
is a initial density increase due to the abrupt temperature 
drop caused by dissipation, and at higher positions, den- 
sity decreases again due to gravity [23] . Density exhibits 
a maximum at z ~ IScr where the density n ~ 0.5, so 
the system cannot be considered as dilute. The density 
profile of the more inelastic particles, B, is plotted as a 
dotted line in Fig la. Its maximum is located at smaller 
z, indicating that they are closer to the bottom of the 
container as compared to the more elastic ones. There- 
fore, particles segregate although the segregation is not 
complete. 

The temperature profiles are also highly inhomoge- 
neous, as shown in Fig. [5^. For both species, tem- 
perature presents a initial abrupt drop, but later (af- 
ter z ~ 20(T both proffies present a linear increase with 
height. This phenomenon was already observed in one- 
component systems, and it is associated to the energy 
transport term, — /zVn, associated to density gradients. 



that appear in granular fluids [23*]. Let us note that the 
maximum density does not coincide with the temperature 
minimum: there is a shift between these two quantities 
which is qualitatively described by a hydrostatic balance 
in presence of gravity .15. , ,23j . 




FIG. 1: Density profiles for two species, density of A, ua, 
(solid line) and the rescaled density of B, ubNa/Nb, (dotted 
line), (a): both species are inelastic with Na = 670, Nb ~ 10, 
OA ~ 0.98, as = 0.7. (b): A is elastic while B is inelastic and 
Na = 640, Nb = 40, and as = 0.5 
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FIG. 2: Temperature profiles for two species, Ta (solid line) 
and Tb (dotted line), (a) and (b) graphics have the same 
parameters as in Fig. [1] 

The described segregation of species A and B is pro- 
duced by their different restitution coefficients as all the 
other properties are the same. To study in more detail 
the effect of the difference of inelasticities and to under- 
stand the origin of this particular segregation, we pro- 
ceed to study the limiting case in which the A particles 
are elastic (a a = 1) and only the B particles are inelastic 
(and consequently, collisions A-B are also inelastic). In 
this way we also limit the parameter space, allowing to a 
more detailed quantitative study. 

Figures [T|3 and show the density and temperature 
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profiles for such case, where particles of type A are elastic 
(q^ — 1) particles B are inelastic (a^ = 0.5) considering 
Nb = 40. It is observed that the main properties of 
the profiles are preserved, even the positive slope of Ta 
despite the A-A collisions are elastic. Partial segregation 
is again observed, where inelastic particles, B, sink to 
the bottom of the container while elastic ones. A, are 
majority at upper layers of the fluid. 

These results are not surprising in view of the predic- 
tions of Ref. [l5| , where it is argued that the segregation 
is produced when the particles with different restitution 
coefficients are immersed in a temperature gradient. The 
gradient is induced by the inelastic collisions, so such 
gradient can be created vibrating a mixture of elastic 
and inelastic particles. The latter ones dissipate some 
of the energy injected by vibration creating a stationary 
state. The hydrodynamic description of the mixture also 
contains the dissipative flux —fiVn, and therefore it is 
expected that the hydrodynamic profiles of density and 
temperature will be equivalent to a full inelastic system. 

Figure [3] shows the temperature ratios Tb{z)/Ta{z) 
for the simulations described in F ig. [21 At low densities 
it was found experimentally (itI . [isj and by employing 
kinetic theory [I6j that such temperature ratio must be 
constant. However, we find a non constant ratio in the z 
direction that only agrees with the result of jT6| at high 
z. Their prediction is valid at low densities in the tracer 
limit, conditions that are only reached in our case for 
large values of z. In the second case, where A par ticles 
are elastic, equivalent predictions were given in [24l. |25| . 
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FIG. 3: Temperature ratios Tb{z)/Ta{z) for the simulations 
described in Fig. [2] Simulation results (solid circles) and 
theoretical prediction in (solid line) for Na = 670, Nb = 
10, OA ~ 0.98, as = 0.7. Simulation results (open squares) 
and theoretical prediction in (dashed line) for A elastic, 
Na = 640, Nb = 40, and as = 0.5 



In order to quantify the segregation, a series of simula- 
tions are performed with Nb ranging from 10 to 160, and 
aB between 0.2 and 0.9. Larger values of Nb or smaller 
restitution coefficients lead to clustering as described in 
[26| . For each simulation we compute the segregation 



parameter, defined as: 

5 — 1 — J dz nA{z)nB{z) j ^ \ dz n\{z) J dz n'^{z) (1) 

where the UAiz) and ub^z) are the local density, as plot- 
ted in Fig. [TJ The segregation parameter is bounded 
between and 1. The value 5—1 corresponds to com- 
plete segregation, as 5 equals 1 only if tia^z) and ub^z) 
do not overlap. On the contrary, 5 = means complete 
mixing, as this value can only be obtained if nB{z) is 
proportional to nA{z). 

The results for 5 are collected in Fig. 0] where the 
quantity 5 is plotted versus the coefficient ub for differ- 
ent values of Nb ■ The fact that 5 is always non vanishing 
confirms that the segregation exists whenever the resti- 
tution coefficients are different. Only in the case when 
as — > 1 the quantity 5 approaches 0, limit in which there 
is no segregation. Note that, for each as, 5 increases 
with Nb- The results confirm that segregation is not 
complete as 5 never gets close to 1. In addition, for each 
simulation, the center of mass of the A and B species are 
computed, Za/Bi finding that Za > Zb in all cases. 




FIG. 4: Segregation parameter 5 as a function of as. Differ- 
ent curves correspond to various concentrations of B particles. 
In the elastic limit a 1 all curves coincide at 5 = 0. 

It can be asked whether the observed segregation could 
compensate the buoyancy force experienced by lighter 
B particles. To verify this idea, a series of simulations, 
keeping fixed aA — 1, o-b = 0.5, Na — 640, Nb = 40, 
and mA{AijjY = f , but varying iriB/rnA is performed. In 
each simulation, the position of the center of masses of 
the A and B species are computed, Za, Zb- The results 
are presented in Fig. [51 The inelastic particles have 
lower center of mass if tob/to^ > 0.37, and therefore 
sinking due to dissipation wins to the buoyancy force in 
this range. On the contrary, buoyancy force dominates 
if iriB/mA < 0.37 and inverse segregation is obtained 
when B particles are lighter than this threshold. The 
segregation parameter 6 does not vanish for any value 
of mB/rriA, indicating that there is no complete mixing 
even at the value of mB/mA = 0.37, where both center 



4 



of masses coincide. The value of 5, however, is minimum 
at this precise mass ratio. 



Z 15 




FIG. 5: Position of the center of mass of A (open squares) and 
B particles (solid circles) as a function of the relative mass 
of B particles rriB/mA- Simulation parameters are fixed to 
OA = 1, as = 0.5, Na = 640, Nb = 40, and mA{Aujf = 1. 
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FIG. 6: Snapshot of the simulated system for Nb = 40, as = 
0.5. A particles are white and B particles black. The border 
of some particles exceed the lateral system size due to the 
periodic boundary conditions in this direction. 



IV. MICROSTRUCTURE 

A snapshot of the previously studied case, aA = 
I, as = 0.5 and Nb = 40 is shown in Fig. [Sj The 
particles of type A are plotted as open circles, and B are 
black symbols. This snapshot suggest that, besides the 
macroscopic partial segregation characterized by i5, there 
is also a micro-segregation, where B particles tend to be 
close to other B particles, differently as if the particles 
where labeled A or B at random. To decide if this ob- 
servation is indeed true, we analyze quantitatively the 
system computing the pair correlation functions. 

As the system is inhomogeneous in the z direction, the 
pair correlation functions became z-dependent. We chose 
to divide the system in horizontal slabs of height lOcr. For 
every pair of particles, the center of mass is obtained and, 
according to it, their relative distance is considered for 
the histogram of distances associated to that particular 
slab. Furthermore, the pairs are classified according to 
the type of particles involved. Finally, the histograms are 
normalized as usual with the local density of each species 
in the slab, such that at long distances the obtained pair 
correlation function approaches to unity. This results in 
functions gfii,{r), for the species /i, = {A,B}, and the 
slab s. 

In most of the studied cases, the region of density 
closest to homogeneous corresponds to the second slab 
(lOtT < z < 20a). In what follows we will present re- 
sults only for this slab and we will suppress accordingly 
the superscript in the pair correlation functions. The 
qualitative properties for other slabs are similar to this 
particular slab, although the effects are reduced because 
at higher slabs the densities are smaller. 



Figure [7] shows the density correlation function for the 
pairs AA and BB for Na = 580, Nb = 100, as = 0.8, 
conditions in which the system develops an average den- 
sity of 71 ~ 0.56. Correlation of distinct particles AB 
have intermediate properties between A A and BB. The 
main noticeable feature of these correlation functions is 
that the first and second peak of the BB function is much 
larger than the A A one. In other words, the large values 
oi gsB at contact means that the are more B-B pairs in 
the system that in a configuration where A and B parti- 
cles are labeled at random. This excess number of B-B 
pairs could be guessed from the figure [HI where one can 
easily locate 4 B-B pairs. 

The structure of the correlation functions in our in- 
elastic system could resemble the radial distribution of 
an elastic fluid at a selected higher density. We tried to 
exploit this idea by comparing the granular distribution 
function with an elastic one by choosing an appropriate 
density. The way to select the density is by adjusting 
the value of the pair distribution function at contact. We 
take this criterion inspired by the kinetic theory of dense 
fluids, where the pair distribution function at contact is 
used to improve the Boltzmann equation including cer- 
tain correlations. With such procedure we flnd different 
fltting densities for AA and BB pairs, being their values: 
nBB = 0.823 and uaa = 0.636. Surprisingly, none of 
them (not even the elastic) agree with the average den- 
sity in the system, n ~ 0.56. This means that the dissi- 
pative particles, type B, are able to modify the structural 
properties of the A particles, which we have chosen to be 
elastic. The fitted elastic pair correlation functions are 
plotted as dashed lines in Fig. [71 Despite of the fitting 
procedure, large discrepancies are observed between our 



5 



correlation functions and the elastic ones. Other fitting 
procedures could have been selected, but none of them 
produce an agreement of the correlation functions in all 
its range. 

The inelastic correlation function shows a strong en- 
hancement of the first peak, followed by a second, and 
at most a third peak, reflecting a short range structure 
of few particle diameters. However, the location of the 
peaks differs form the elastic correlation function, indi- 
cating changes in the microstructure. The fast decay of 
the correlation function in inelastic systems has been ob- 
served in homogeneously heated inelastic gases as well 
[131 . This analysis confirms the existence of a microseg- 
regation of particles of type B. 





FIG. 7: Solid lines: radial distribution functions gAA{r), and 
gBBif). Nb = 100, as = 0.8. Dashed lines: equilibrium 
radial distribution function given in Ref . [l^l , with the density 
fitted to give the same value of qaa or qbb at contact (left 
UAA = 0.636, right ubb = 0.823) 

In addition, the pair correlation function at contact 
X = 5(0'^) is computed for the pairs AA, AB, and BB. 
The values of x grow when increasing Nb and decreasing 
a_B, obeying that xbb > Xab > Xaa- Finally let us re- 
mark that Xbb can reach very high values. For instance, 
at Nb = 160 and ub = 0.5, xs — 60, despite that the 
average total density is only n ~ 0.5. 

It has been proposed [2^, based on an Enskog-like ki- 
netic model that, due to the smaller outgoing velocity 
after a collision, particles stay longer in the close vicinity, 
when the restitution coefficient is less than one. These 
arguments leads to an enhancement of the pair correla- 
tion function in terms of the restitution coefficient a, as 



(2) 



Here xo is the pair correlation function of elastic parti- 
cles at the same density, that are given by the Verlet- 
Levesque [S^] and Carnahan-Starling [3l| factors for 2 
and 3 dimensions. 

To test if this prediction is valid we have plotted in Fig. 
[5] the ratios xbb/xaa against as and xab/ Xaa against 



aAB and compared them with the factors {I+ub) / {'i.aB) 
and {1 + uab) I C^olab) respectively. Dividing by xaa we 
get rid of the xo factor in The figure indicates that 
the Enskog model is not enough to describe the high val- 
ues of the pair correlation when one inelastic particle is 
present, as was already quoted by its author [2^ . A possi- 
ble origin of this effect are recoUisions, that are not taken 
into account in Enskog's theory, a phenomenon that is 
pronounced in granular systems. Inelastic particles, af- 
ter a collision, separate at a slower speed, increasing the 
possibility of having a collision with a third particle, ap- 
proaching it again the original pair. In some way, inelas- 
ticity increases the so-called cage effect in liquids. 
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FIG. 8: Contact pair correlation function quotients as a func- 
tion of the restitution coefficient, xbb /xaa results from sim- 
ulations (filled circles) and the Enskog theoretical value (solid 
line); and xab/ Xaa results from simulations (open trian- 
gles) and the Enskog theoretical value (dotted line). Here 
Nb = 120. 

Besides density correlations induced by B particles, 
they also create a local decrease of temperature because 
their collisions are inelastic. We have computed the ra- 
dial temperature of particles of type ^ located at a dis- 
tance r oi & V particle: T^^{r). Figure [3] shows these 
radial temperature functions. Their asymptotic values 
are the average temperatures of A and B particles in the 
slab 10a < r < 20cr. Besides that B particles are colder 
than A particles, they produce a local decrease of tem- 
perature both for A and B particles. This effect is more 
pronounced when increasing the inelasticity but almost 
independent of the concentration of B particles, whose 
main effect is to reduce globally the temperature. As 
before this is a local effect that extends for about 2cr. 



V. ORIGIN OF THE GLOBAL SEGREGATION 

So far we have argued along the paper that B particles 
locally modify the structure of B particles (inelastic) but 
also of A particles (elastic). For instance, Fig. [7] shows 
that the density around B particles is higher than around 
A particles. The excess mass around B particles with 
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FIG. 9: Radial temperatures T^i^ of a particle of type ^ 
around a particle of type v. Nb = 40, qs = 0.5. 

respect A particles is defined as, 

5in= dr2TTr{pB{r)-pA{r)), (3) 

J a 

where (r) is the total density around a particle of type 
[I. The excess mass is positive for all values of Na^Nb 
and restitution coefficient ub- Similarly, the local tem- 
perature around a B particle is lower than the local tem- 
perature around a A particle (see Fig. [9|). 

These two results imply that around a B particle it is 
developed a dense and cold region, that we call a cold 
droplet. The characteristic size of the droplet is micro- 
scopic; in agreement with the plots of Sect IV, must be 
about 2-3 a. The development of a cold droplet around 
each B particle increases their effective mass and also, due 
to the lower temperature, decreases the pressure. This 
mechanism could be similar (at the local level) to the 
clustering instability described in [s^]. In our case the 
unstable process does not continue as the energy injec- 
tion due to the vibration breaks the droplet, that must 
be seen as dynamical object, that forms and evanesces 
continuously. Once the cold droplet is formed around 
B the buoyancy force in weaker than its effective weight 
and therefore there is a net force that tends to sink the 
cluster and consequently the B particle. 

This continuously sink of the B particle is the respon- 
sible of the macroscopic segregation described in Sect. 
III. 



VI. SUMMARY AND CONCLUSIONS 

The main conclusion of the present paper is that dif- 
ferent restitution coefficients alone create segregation in 
a binary mixture vertically vibrated. The restitution co- 
efficients must be considered, in addition to the usual 
material properties (mass ratio and diameter ratio) in 
order to describe accurately the segregation. The effect 



of the inelasticity is such that the most inelastic particles 
sink to the bottom of the container while the less inelas- 
tic ones raise to the top. Segregation is not complete, 
however, but only partial. The density profiles of each 
species shows a maximum which is located in a differ- 
ent position depending on the inelasticities. Concerning 
the temperatures, most dissipative particles have a lower 
temperature than the most elastic ones, both having the 
characteristic shape of vibrofluidized system (fast cooling 
away from the moving boundary followed by a heating 
that grows linearly with the distance). 

The segregation effects presented here also appear by 
vibrating a mixture of elastic and inelastic particles. 
Again inelastic particles migrate to the bottom of the 
container and elastic one prefer the upper part. The 
temperature distribution also looks similar to the full in- 
elastic mixture. 

Besides the macroscopic segregation there is also seg- 
regation at the microscopic scale. With the word micro- 
scopic we refer to properties at distances of few particles 
diameters. In our case we find a notorious increase of the 
probability of finding two inelastic particles together as 
compared with the less inelastic or elastic ones. This en- 
hancement cannot be described by only considering kine- 
matic properties. Our guess is that dynamic correlations 
are required to properly describe these correlations. 

Finally, we propose a mesoscopic explanation to the 
segregation: the most inelastic particles induce locally a 
region of high density and low temperature, resembling a 
cold droplet that falls in a gravitational field. The droplet 
is created by the dissipation in a way that resembles the 
clustering instability of the granular gases. 

The effect of the different restitution coefficients may 
act in an opposite direction as the usual Brazil nut effect. 
For instance, consider a vibrated granular fluid and insert 
a large (or light) intruder than tends to move upward. If 
the intruder is very dissipative, it will move downwards, 
as we have described in this paper. Which is the final 
position? Which force does finally win? Partial answer 
is presented in Fig. [5l where we show the competition 
between the buoyancy force and the sinking effect due 
to dissipation. Moreover, could two large an inelastic 
intruders come together by the effect of inelasticity alone, 
as shown in the Sect. IV of the present paper? Further 
research is needed in order to answer these questions. 
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